#install.packages(c("lidR","terra", "dplyr", "viridis))
library(lidR)
library(viridis)
library(dplyr)
library(terra)uav_lidar_inventory pt. 3
LiDAR Forest Inventory - Part 3
Main contents of this exercise
Compare Terrestrial Laser Scanning (TLS)-based and Unmanned Aerial Vehicle (UAV)-based point clouds for forest inventries.
Compare UAV and TLS-based instance segmentations to a full forest inventory of the ECOSENSE site (Ettenheim, Germany).
Instance segmentations delineate individual trees in point clouds. Here, we test SegmentAnyTree. A small visual recap on SegmentAnyTree, a tool to segment individual trees in LiDAR point clouds across different LiDAR data types:
Moreover, we will look into AI-based tree species recognition that was applied on the invidiual tree segmentations.
Remember that we looked into how the point clouds of ULS and TLS differ from one another, here is a quick reminder which shows the two side by side:
Load Packages, paths, etc…
Reminder: this script assumes that your Quarto document is in the same directory as your data. If this is not the case, consider to set your working directory with setwd("path/to/your/data/folder/").
Load & inspect the data
Inventory data
Load the tree inventory data from ECOSENSE (full inventory).
*.gpkg stands for geopackage format. We can load that in R using terra as SpatVector object.
tree_inv <- vect("tree_inventory_final2.gpkg")head(tree_inv) species tls_name tls_TH tls_DBH tls_X tls_Y
1 BE T2792 27.28962149 0.088146948 416700.687 5346611.676
2 BE T2756 27.65825122 0.16474044 416698.162 5346611.911
3 BE T2763 26.53824404 0.136313621 416699.4069 5346613.154
4 BE T2793 28.58093068 0.363458403 416704.5793 5346614.498
5 BE T2816 28.73643552 0.387836247 416703.9667 5346609.953
6 BE <NA> <NA> <NA> <NA>
tls_Longitude tls_Latitude odk_add_tree.loop_time
1 7.877502550380541 48.26723262917562 2024-11-26T11:23:27Z
2 7.877468485304306 48.26723441090238 2024-11-26T11:36:16Z
3 7.877485012699943 48.26724575547465 2024-11-26T11:41:34Z
4 7.877554434593113 48.26725852509873 2024-11-26T11:43:50Z
5 7.877547076487936 48.26721756208074 2024-11-26T11:46:38Z
6 <NA> <NA> 2024-11-26T11:48:11Z
odk_add_tree.plot_id odk_add_tree.tree_id odk_add_tree.new_tree
1 <NA> <NA> true
2 16 <NA> true
3 <NA> <NA> true
4 <NA> <NA> true
5 <NA> <NA> true
6 <NA> <NA> true
odk_qrcode.qrcode_id odk_qrcode.qrcode_bu_id
1 https://dt.unr.uni-freiburg.de/ecosense/16_27
2 05865987 16_28
3 https://dt.unr.uni-freiburg.de/ecosense/16_29
4 https://dt.unr.uni-freiburg.de/ecosense/16_30
5 https://dt.unr.uni-freiburg.de/ecosense/16_31
6 https://dt.unr.uni-freiburg.de/ecosense/16_32
odk_new_tree_details.circ odk_new_tree_details.other_species
1 28 <NA>
2 52 <NA>
3 41 <NA>
4 115 <NA>
5 125 <NA>
6 25 <NA>
odk_tree_img odk_gps_tree.Latitude odk_gps_tree.Longitude
1 1732620954736.jpg 48.2672236483333 7.877503075
2 1732621157835.jpg 48.2672254283333 7.87747073833333
3 1732621415109.jpg 48.2672363233333 7.87748745333333
4 1732621558170.jpg 48.2672479383333 7.87755160333333
5 1732621683746.jpg 48.2672075166667 7.87754983833333
6 1732621738341.jpg 48.2672061566667 7.87752802833333
odk_gps_tree.Altitude odk_gps_tree.Accuracy odk_tree_comment
1 526.659 0.01
2 526.557 0.01 No barcode tag
3 526.663 0.01
4 526.798 0.01
5 526.669 0.01
6 526.925 0.01
odk_PARENT_KEY
1 uuid:e209a154-c01c-4534-b144-cfa7d7713b22
2 uuid:e209a154-c01c-4534-b144-cfa7d7713b22
3 uuid:e209a154-c01c-4534-b144-cfa7d7713b22
4 uuid:e209a154-c01c-4534-b144-cfa7d7713b22
5 uuid:e209a154-c01c-4534-b144-cfa7d7713b22
6 uuid:e209a154-c01c-4534-b144-cfa7d7713b22
odk_KEY odk_diameter
1 uuid:e209a154-c01c-4534-b144-cfa7d7713b22/loop[1] 0.0891267681314614
2 uuid:e209a154-c01c-4534-b144-cfa7d7713b22/loop[2] 0.165521140815571
3 uuid:e209a154-c01c-4534-b144-cfa7d7713b22/loop[3] 0.130507053335354
4 uuid:e209a154-c01c-4534-b144-cfa7d7713b22/loop[4] 0.366056369111359
5 uuid:e209a154-c01c-4534-b144-cfa7d7713b22/loop[5] 0.397887357729738
6 uuid:e209a154-c01c-4534-b144-cfa7d7713b22/loop[6] 0.0795774715459477
odk_shift_distance_m odk_shifted_longitude odk_shifted_latitude
1 1.08912676813146 7.877502 48.26723
2 1.16552114081557 7.877468 48.26723
3 1.13050705333535 7.877485 48.26725
4 1.36605636911136 7.877555 48.26726
5 1.39788735772974 7.877547 48.26722
6 1.07957747154595 7.877526 48.26722
dbh_difference new_longitude new_latitude
1 9.798201314614036e-4 7.877502550380541 48.26723262917562
2 7.807008155710227e-4 7.877468485304306 48.26723441090238
3 0.00580656766464599 7.877485012699943 48.26724575547465
4 0.00259796611135904 7.877554434593113 48.26725852509873
5 0.010051110729737966 7.877547076487936 48.26721756208074
6 <NA> 7.87752802833322 48.2672159464496
geometry last_trees TreeID_test_area
1 c(416700.727386307, 5346611.77496821) 16_27
2 c(416698.331770989, 5346612.08489961) 05865987
3 c(416699.589345929, 5346613.24243561) 16_29
4 c(416704.372125069, 5346614.70123457) 16_30
5 c(416704.175926162, 5346610.24241929) 16_31
6 c(416702.550551442, 5346609.7940927) 16_32
layer path TreeID
1 last_trees_1 — last_trees last_trees_1.gpkg|layername=last_trees 16_27
2 last_trees_1 — last_trees last_trees_1.gpkg|layername=last_trees 16_28
3 last_trees_1 — last_trees last_trees_1.gpkg|layername=last_trees 16_29
4 last_trees_1 — last_trees last_trees_1.gpkg|layername=last_trees 16_30
5 last_trees_1 — last_trees last_trees_1.gpkg|layername=last_trees 16_31
6 last_trees_1 — last_trees last_trees_1.gpkg|layername=last_trees 16_32
species_name
1 Fagus_sylvatica
2 Fagus_sylvatica
3 Fagus_sylvatica
4 Fagus_sylvatica
5 Fagus_sylvatica
6 Fagus_sylvatica….puh, many attributes. For sure not all of them are useful for this analysis. We will mostly focus on the species, the DBH and of course the geolocation:
plot(tree_inv, col = "red")
# please ignore
#tree_inv_spec <- read.csv("lookup_odk2.csv", sep = ",", header = T)
#tree_inv2 <- merge(tree_inv, tree_inv_spec, by.x = "species", by.y = "abbr")
#writeVector(tree_inv2, "D:/output_file.gpkg", filetype = "GPKG")TLS & UAV-based inventory data
UAV-based point clouds with AI-based tree segmentation
pc_uav <- readLAS("uav_pointcloud_red.las")
pc_uavclass : LAS (v1.4 format 6)
memory : 326.8 Mb
extent : 416648.3, 416838.9, 5346560, 5346871 (xmin, xmax, ymin, ymax)
coord. ref. : NA
area : 41168 units²
points : 6.12 million points
density : 148.66 points/units²TLS-based point clouds with AI-based tree segmentation
pc_tls <- readLAS("tls_pointcloud_red.las")
pc_tlsclass : LAS (v1.4 format 6)
memory : 294.8 Mb
extent : 416648.2, 416839, 5346560, 5346871 (xmin, xmax, ymin, ymax)
coord. ref. : NA
area : 41140 units²
points : 5.52 million points
density : 134.17 points/units²Let´s plot the data for a visual comparison:
col <- sample(viridis(2070, option = "mako"))
p <- plot(pc_uav, color = "PredInstance", pal = col, legend = F, axis = T, nbreaks = length(col))col <- sample(viridis(2070, option = "mako"))
p <- plot(pc_tls, color = "PredInstance", pal = col, legend = F, axis = T, nbreaks = length(col))UAV-based segmentation vs full inventory
Remember how we extracted the tree cordinates in the previous exercise? We will apply this again.
Normalize point cloud on Z-axis
mycsf <- csf(TRUE, 1, 1, time_step = 1)
pc_uav <- classify_ground(pc_uav, mycsf)
pc_uav_n <- normalize_height(pc_uav, knnidw())Inverse distance weighting: [=====================================-------------] 74% (8 threads)
Inverse distance weighting: [=====================================-------------] 75% (8 threads)
Inverse distance weighting: [======================================------------] 76% (8 threads)
Inverse distance weighting: [======================================------------] 77% (8 threads)
Inverse distance weighting: [=======================================-----------] 78% (8 threads)
Inverse distance weighting: [=======================================-----------] 79% (8 threads)
Inverse distance weighting: [========================================----------] 80% (8 threads)
Inverse distance weighting: [========================================----------] 81% (8 threads)
Inverse distance weighting: [=========================================---------] 82% (8 threads)
Inverse distance weighting: [=========================================---------] 83% (8 threads)
Inverse distance weighting: [==========================================--------] 84% (8 threads)
Inverse distance weighting: [==========================================--------] 85% (8 threads)
Inverse distance weighting: [===========================================-------] 86% (8 threads)
Inverse distance weighting: [===========================================-------] 87% (8 threads)
Inverse distance weighting: [============================================------] 88% (8 threads)
Inverse distance weighting: [============================================------] 89% (8 threads)
Inverse distance weighting: [=============================================-----] 90% (8 threads)
Inverse distance weighting: [=============================================-----] 91% (8 threads)
Inverse distance weighting: [==============================================----] 92% (8 threads)
Inverse distance weighting: [==============================================----] 93% (8 threads)
Inverse distance weighting: [===============================================---] 94% (8 threads)
Inverse distance weighting: [===============================================---] 95% (8 threads)
Inverse distance weighting: [================================================--] 96% (8 threads)
Inverse distance weighting: [================================================--] 97% (8 threads)
Inverse distance weighting: [=================================================-] 98% (8 threads)
Inverse distance weighting: [=================================================-] 99% (8 threads)
Inverse distance weighting: [==================================================] 100% (8 threads)plot(pc_uav_n, axis = T)Cut a horizonal transect for stem extraction (we slice the point cloud close to the ground).
pc_uav_n_trans <- filter_poi(pc_uav_n, Z >= 0.1 & Z <= 4)
col <- sample(viridis(2070, option = "mako"))plot(pc_uav_n_trans, color = "PredInstance", pal = col, legend = F, axis = T, nbreaks = length(col))Derive the centroids for stems (as proxy for the tree position):
# Calculate centroids for each instance
centroids_uav <- pc_uav_n_trans@data %>%
group_by(PredInstance) %>%
summarize(
X = mean(X),
Y = mean(Y),
Z = 0,
species = as.integer(names(which.max(table(species_id)))),
.groups = "drop" # Avoid unnecessary grouping in the result
)
# Print centroids
print(centroids_uav)# A tibble: 1,751 × 5
PredInstance X Y Z species
<int> <dbl> <dbl> <dbl> <int>
1 0 416745. 5346740. 0 0
2 1 416779. 5346568. 0 19
3 2 416772. 5346570. 0 10
4 5 416784. 5346569. 0 10
5 6 416784. 5346571. 0 10
6 7 416780. 5346569. 0 10
7 8 416787. 5346566. 0 25
8 9 416796. 5346561. 0 19
9 10 416790. 5346563. 0 10
10 12 416797. 5346568. 0 10
# ℹ 1,741 more rowsConvert the centroids to a geospatial vector format (SpatVector), so we can compare it with the inventory data (which is also in vector format):
centroids_uav_sv <- vect(centroids_uav, geom = c("X", "Y"), crs = "EPSG:32632")Visually compare the UAV-based centroids with the inventory data
plot(centroids_uav_sv)
points(tree_inv, col = "red")
Quantitatively compare the UAV-based centroids with the inventory data
buffer_tree_inv <- buffer(tree_inv, width = 1)
plot(buffer_tree_inv)
points(tree_inv, col = "red")
matches_uav <- intersect(buffer_tree_inv, centroids_uav_sv)plot(tree_inv, col = "red")
points(matches_uav, col = "blue", cex = 1.3)
nrow(matches_uav)/nrow(tree_inv)*100[1] 66.16786Is there a systematic bias in the detection?
unmatched_uav <- buffer_tree_inv[!(buffer_tree_inv$TreeID %in% matches_uav$TreeID), ]par(mfrow = c(1,2))
hist(as.numeric(matches_uav$tls_DBH),
xlim = c(0,0.9),
ylim=c(0,200),
main = "DBH matched",
xlab = "DBH [cm]")
hist(as.numeric(unmatched_uav$tls_DBH),
xlim = c(0,0.9),
ylim=c(0,200),
main = "DBH unmatched",
xlab = "DBH [cm]")
Not all trees are matched, but at least there seems to be no strong bias for small/large trees.
TLS-based segmentation vs full inventory
Same procedure as for the UAV data. First, we normalize the point cloud…
mycsf <- csf(TRUE, 1, 1, time_step = 1)
pc_tls <- classify_ground(pc_tls, mycsf)
pc_tls_n <- normalize_height(pc_tls, knnidw())Inverse distance weighting: [=======================================-----------] 79% (8 threads)
Inverse distance weighting: [========================================----------] 80% (8 threads)
Inverse distance weighting: [========================================----------] 81% (8 threads)
Inverse distance weighting: [=========================================---------] 82% (8 threads)
Inverse distance weighting: [=========================================---------] 83% (8 threads)
Inverse distance weighting: [==========================================--------] 84% (8 threads)
Inverse distance weighting: [==========================================--------] 85% (8 threads)
Inverse distance weighting: [===========================================-------] 86% (8 threads)
Inverse distance weighting: [===========================================-------] 87% (8 threads)
Inverse distance weighting: [============================================------] 88% (8 threads)
Inverse distance weighting: [============================================------] 89% (8 threads)
Inverse distance weighting: [=============================================-----] 90% (8 threads)
Inverse distance weighting: [=============================================-----] 91% (8 threads)
Inverse distance weighting: [==============================================----] 92% (8 threads)
Inverse distance weighting: [==============================================----] 93% (8 threads)
Inverse distance weighting: [===============================================---] 94% (8 threads)
Inverse distance weighting: [===============================================---] 95% (8 threads)
Inverse distance weighting: [================================================--] 96% (8 threads)
Inverse distance weighting: [================================================--] 97% (8 threads)
Inverse distance weighting: [=================================================-] 98% (8 threads)
Inverse distance weighting: [=================================================-] 99% (8 threads)
Inverse distance weighting: [==================================================] 100% (8 threads)plot(pc_tls_n, axis = T)… create a transect of stems…
pc_tls_n_trans <- filter_poi(pc_tls_n, Z >= 0.1 & Z <= 4)
col <- sample(viridis(2070, option = "mako"))plot(pc_tls_n_trans, color = "PredInstance", pal = col, legend = F, axis = T, nbreaks = length(col))… find the centroids of the stems.
# Calculate centroids for each instance
centroids_tls <- pc_tls_n_trans@data %>%
group_by(PredInstance) %>%
summarize(
X = mean(X),
Y = mean(Y),
Z = 0,
species = as.integer(names(which.max(table(species_id)))),
.groups = "drop" # Avoid unnecessary grouping in the result
)
# Print centroids
print(centroids_tls)# A tibble: 2,093 × 5
PredInstance X Y Z species
<int> <dbl> <dbl> <dbl> <int>
1 0 416740. 5346724. 0 0
2 797 416718. 5346589. 0 10
3 798 416723. 5346588. 0 10
4 802 416727. 5346587. 0 10
5 803 416735. 5346586. 0 10
6 807 416736. 5346582. 0 10
7 808 416730. 5346586. 0 10
8 809 416737. 5346583. 0 3
9 810 416734. 5346583. 0 10
10 811 416742. 5346580. 0 3
# ℹ 2,083 more rowsLet´s compare the result to the inventory data:
centroids_tls_sv <- vect(centroids_tls, geom = c("X", "Y"), crs = "EPSG:32632")plot(centroids_tls_sv)
points(tree_inv, col = "red")
buffer_tree_inv <- buffer(tree_inv, width = 1)
plot(buffer_tree_inv)
points(tree_inv, col = "red")
matches_tls <- intersect(centroids_tls_sv, buffer_tree_inv)plot(buffer_tree_inv)
points(matches_tls, col = "blue")
nrow(matches_tls)/nrow(tree_inv)*100[1] 81.52245TLS & UAV tree height comparison
Above we have looked at the coordinates. But we know that we can extract much more, such as the height.
Extract maximum Z value per instance for both TLS and UAV data:
# Find the highest point for each instance
uav_tree_heights <- pc_uav_n@data %>%
group_by(PredInstance) %>%
slice_max(Z, n = 1) %>% # Get the row with the maximum Z value for each instance
ungroup()
uav_tree_heights <- uav_tree_heights[, c("X", "Y", "Z", "PredInstance")]
# View the results
print(uav_tree_heights)# A tibble: 2,161 × 4
X Y Z PredInstance
<dbl> <dbl> <dbl> <int>
1 416754. 5346752. 46.6 0
2 416778. 5346567. 28.4 1
3 416772. 5346569. 24.7 2
4 416778. 5346569. 28.6 4
5 416788. 5346570 13.1 5
6 416788. 5346571. 30.4 6
7 416780. 5346570. 29.2 7
8 416788. 5346566. 34.2 8
9 416796. 5346561. 27.0 9
10 416791. 5346562. 26.7 10
# ℹ 2,151 more rows# Find the highest point for each instance
tls_tree_heights <- pc_tls_n@data %>%
group_by(PredInstance) %>%
slice_max(Z, n = 1) %>% # Get the row with the maximum Z value for each instance
ungroup()
tls_tree_heights <- tls_tree_heights[, c("X", "Y", "Z", "PredInstance")]
# View the results
print(tls_tree_heights)# A tibble: 2,578 × 4
X Y Z PredInstance
<dbl> <dbl> <dbl> <int>
1 416750. 5346755. 37.9 0
2 416711. 5346592. 15.7 669
3 416760. 5346574. 19.0 700
4 416756. 5346575. 18.3 703
5 416768. 5346570. 24.4 709
6 416785. 5346564. 27.6 717
7 416699. 5346597. 11.8 787
8 416692. 5346599. 18.6 789
9 416716. 5346590. 28.7 797
10 416719. 5346589. 26.4 798
# ℹ 2,568 more rowsLet´s compare this visually using histograms:
par(mfrow = c(1,2))
hist(tls_tree_heights$Z, ylim = c(0,900), main = "TLS-based height distribution")
hist(uav_tree_heights$Z, ylim = c(0,900), main = "UAV-based height distribution")
paste("number TLS-based trees: ", nrow(tls_tree_heights))[1] "number TLS-based trees: 2578"print("TLS-based tree heights")[1] "TLS-based tree heights"summary(tls_tree_heights$Z) Min. 1st Qu. Median Mean 3rd Qu. Max.
0.00 10.56 23.44 19.51 27.96 46.18 paste("number UAV-based trees: ", nrow(uav_tree_heights))[1] "number UAV-based trees: 2161"print("UAV-based tree heights")[1] "UAV-based tree heights"summary(uav_tree_heights$Z) Min. 1st Qu. Median Mean 3rd Qu. Max.
2.065 12.636 25.462 21.584 28.962 48.443 When looking at the transects … can you explain the results?
AI-based Species identification
We will test
DetailViewfrom Julian Frey et al.: https://github.com/JulFrey/DetailViewYou can read on the details here: Puliti, S., Lines, E. R., Müllerová, J., Frey, J., Schindler, Z., Straker, A., ... & Astrup, R. (2024). Benchmarking tree species classification from proximally-sensed laser scanning data: introducing the FOR-species20K dataset. arXiv preprint arXiv:2408.06507.
Detail view requires point clouds of invidiual trees. Thus, we are perfectly prepared as with SegmentAnyTree, we already extracted single trees. Here you can see some of the training data of DetailView:
Currently, DetailView is one of the most accurate tree species classifiaction modelds for LiDAR data:
Let´s first have look at the reality (forest inventory data of ECOSENSE):
par(mar = c(5, 15, 2, 5))
barplot(sort(table(tree_inv$species_name)), las = 1, horiz = T)
The UAV- and TLS-based predictions come in species-codes, not species names. So we first load a look up table to attach (merge) the actual names to our data:
species_ai <- read.csv("lookup_ai.csv")… attach the species name to the point clouds….
pc_tls@data <- merge(pc_tls@data, species_ai, by = "species_id", all.x = TRUE)
pc_uav@data <- merge(pc_uav@data, species_ai, by = "species_id", all.x = TRUE)… and to the extract stem coordinates (centroids):
centroids_tls_sv2 <- merge(centroids_tls_sv , species_ai, by.x ="species", by.y = "species_id", all.x = TRUE)
centroids_uav_sv2 <- merge(centroids_uav_sv , species_ai, by.x ="species", by.y = "species_id", all.x = TRUE)Visual comparison of the species classification
p <- plot(pc_tls, color = "species_id", legend = T, axis = T)
p <- plot(pc_uav, color = "species_id", legend = T, axis = T)Okay, the results “look” already quite different. Let´s check the stats:
uav_n_spec <- pc_uav@data %>%
as.data.frame() %>%
group_by(species) %>%
summarise(number_of_instances = n_distinct(PredInstance)) %>%
arrange(desc(number_of_instances))
uav_n_spec# A tibble: 15 × 2
species number_of_instances
<chr> <int>
1 Fagus_sylvatica 1110
2 Betula_pendula 408
3 Pinus_radiata 185
4 Carpinus_betulus 91
5 Pseudotsuga_menziesii 83
6 Quercus_rubra 77
7 Pinus_sylvestris 68
8 Picea_abies 63
9 Acer_pseudoplatanus 34
10 Acer_saccharum 16
11 Abies_alba 12
12 Acer_campestre 6
13 Quercus_petraea 4
14 Larix_decidua 3
15 Quercus_robur 1tls_n_spec <- pc_tls@data %>%
as.data.frame() %>%
group_by(species) %>%
summarise(number_of_instances = n_distinct(PredInstance)) %>%
arrange(desc(number_of_instances))
tls_n_spec# A tibble: 23 × 2
species number_of_instances
<chr> <int>
1 Fagus_sylvatica 1024
2 Acer_saccharum 626
3 Pseudotsuga_menziesii 196
4 Abies_alba 95
5 Picea_abies 92
6 Pinus_sylvestris 56
7 Quercus_robur 52
8 Corylus_avellana 46
9 Betula_pendula 45
10 Carpinus_betulus 43
# ℹ 13 more rowssort(table(pc_tls$species))
par(mar = c(10, 5, 2, 5))
# Convert species to a named vector
species_counts <- setNames(tls_n_spec$number_of_instances, tls_n_spec$species)
# Create the barplot
barplot(species_counts,
las = 2, # Rotate axis labels for readability
col = "skyblue", # Add some color
main = "Species numbers from TLS data", # Add a title
ylab = "Number of Instances") # Add Y-axis label
par(mar = c(10, 5, 2, 5))
# Convert species to a named vector
species_counts <- setNames(uav_n_spec$number_of_instances, uav_n_spec$species)
# Create the barplot
barplot(species_counts,
las = 2, # Rotate axis labels for readability
col = "skyblue", # Add some color
main = "Species numbers from UAV data", # Add a title
ylab = "Number of Instances") # Add Y-axis label
Compare the tree species distribution among datasets:
# Extract unique species across all SpatVectors
all_species <- unique(c(tree_inv$species_name, centroids_tls_sv2$species.y, centroids_uav_sv2$species.y))
tree_inv$species.y <- tree_inv$species_name
# Create a color palette for all species
species_colors <- setNames(rainbow(length(all_species)), all_species)
# Function to assign colors based on species
get_colors <- function(spatvector, color_mapping) {
color_mapping[spatvector$species.y]
}
# Adjust layout to include space for the legend
par(mfrow = c(1, 3), mar = c(4, 4, 2, 1)) # 2x2 grid, adjust margins as needed
plot(tree_inv,
col = get_colors(tree_inv, species_colors),
main = "Inventory",
legend = FALSE)
plot(centroids_tls_sv2,
col = get_colors(centroids_tls_sv2, species_colors),
main = "Species by TLS",
legend = FALSE)
plot(centroids_uav_sv2,
col = get_colors(centroids_uav_sv2, species_colors),
main = "Species by UAV",
legend = FALSE)
par(mfrow = c(1, 3), mar = c(1, 4, 1, 1))
plot.new() # Create an empty plot for the legend
legend("center", legend = names(species_colors), fill = species_colors, title = "Species")
Open Questions & Outlook
How could we potentially improve the tree species classification?
Can we trust UAV-based or TLS-based inventories?
Can we trust field-based inventories in terms of coverage, sampling and measurement uncertainty?
What are related advantages or disadvantages?
How can we move on with UAV- or TLS-based inventories to estimate basal area, biomass or timber volume?
How will forest inventories look like in 2050?